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Taking into account the fact that overload failures in real-world functional networks are usually 
caused by extreme values of temporally fluctuating loads that exceed the allowable range, we study 
the robustness of scale-free networks against cascading overload failures induced by fluctuating loads. 

In our model, loads are described by random walkers moving on a network and a node fails when 
the number of walkers on the node is beyond the node capacity. Our results obtained by using 
the generating function method shows that scale-free networks are more robust against cascading 
overload failures than Erdos-Renyi random graphs with homogeneous degree distributions. This 
conclusion is contrary to that predicted by previous works which neglect the effect of fluctuations 
of loads. 


I. INTRODUCTION 

Modern human societies are supported by various func¬ 
tional networks, such as power grids, the Internet, road 
systems, and corporate transaction networks fli- Since 
the function of a network is guaranteed by its global 
connectivity, a decomposition of the network into dis¬ 
connected components by failures on network elements 
(nodes or links) might induce the breakdown of the func¬ 
tion, which causes a fatal damage to our daily life. It is 
thus crucial to elucidate what type of network topology 
is resilient to failures. In this context, the vulnerability of 
complex networks against random failures and targeted 
attacks has been extensively studied [l-@ . In these stud¬ 
ies, failed nodes are removed from a network at the same 
time. In addition to such simultaneous failures, many 
of actual breakdowns are driven by a chain of failures 
triggered by initial failures of a single or a few nodes. 
Initial failures and subsequent avalanche of failures in 
functional networks are often induced by loads exceeding 
node capacities. For instance, a bankruptcy of an insol¬ 
vent company and subsequent chain bankruptcies due to 
redistributed debts of the bankrupt company can be re¬ 
garded as a process of such cascading overload failures in 
a corporate transaction network. 

In order to describe how networks lose their global con¬ 
nectivity by cascading overload failures, Motter and Lai 
proposed a model in which the load at a node is given by 
the betweenness centrality of the node Q. Their model 
predicts that the initial removal of the highest degree (or 
highest load) node leads a large scale cascade and scale- 
free networks are more fragile against cascading overload 
failures than homogeneous networks. The vulnerability 
of scale-free networks has also been found in other models 
of cascading failures [l0l - [T5j |. In these models, overload 
failures are caused by non-fluctuating loads determined 
uniquely by the network structure that exceed the ca¬ 
pacity of each node. It is, however, quite general that 
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the load on a node fluctuates temporally [l6l - [l8j| . An 
overload failure takes place when an instantaneous value 
of fluctuating load exceeds the capacity, as in the cases 
of river flooding or financial collapse in a country M- 
The property of cascading overload failures induced by 
extreme events might be largely different from that by 
non-fluctuating loads. 

However, the robustness of networks against cascading 
failures induced by fluctuating loads has not been widely 
argued so far [2^ . Among several ways to describe fluctu¬ 
ating loads, Kishore et al. modeled them by random walk¬ 
ers on a network [lil, |22| ■ They calculated the overload 
probability that the number of random walkers on a node 
exceeds the predetermined node capacity. Moreover, ap¬ 
plying this theory, the network robustness against non¬ 
cascading overload failures has been studied, in which 
nodes are simultaneously removed once according to the 
overload probability [^. Although these theories 
take into account the temporal fluctuations of loads, cas¬ 
cade processes triggered by the initial overload failures 
have not been considered. In this paper, we examine the 
robustness of complex networks against cascading over¬ 
load failures induced by the extreme value of fluctuating 
loads. Adopting the random walker model proposed by 
Kishore et al. and their theory of the overload probability 
[U [ 23 , we present a simple model to describe cascades 
of overload failures caused by fluctuating loads and cal¬ 
culate the size of the giant component after completing 
the cascade by using the generating function formalism. 
Our main result shows that scale-free networks are more 
robust against cascading overload failures than homoge¬ 
neous networks, which is contrary to that predicted by 
previous works 

The rest of this paper is organized as follows. In Sec.llll 
we present a model to describe cascading overload fail¬ 
ures on a complex network based on the random walker 
model proposed by Kishore et al. [U |2^ . In Sec. Illli 
we explain the method to calculate the size of the giant 
component after completing the cascade process by uti¬ 
lizing the master equation for the probability of a node 
to have the initial degree ko and the degree A: at a cas¬ 
cade step T and the generating function formalism. Our 
results are presented in Sec. lYl Section El is devoted to 
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the summary and concluding remarks. which Wq random walkers exist, and determine the ca¬ 

pacity qk of each node according to Eq. ©• 

(ii) At each time step r, assign WV random walkers to 
II. MODEL the network Qr at time r. The total load Wr is given by 


In a functional network, some sort of “flow” is often 
required to realize its functionality, and at the same time 
the flow plays a role of a “load” in the network, such as 
electric current in a power grid or packet transfer on the 
Internet. The relation between the average and fluctua¬ 
tion of flow in such networks has been investigated empir¬ 
ically and theoretically dMl. These studies elucidated 
that the flux fluctuations at a node have the same scaling 
behavior with the fluctuations of the number of random 
walkers on the node. Inspired by this fact, Kishore et 
al. modeled fluctuating loads by random walkers moving 
on a network [U, , where the number of walkers in¬ 

dicates the amount of loads. The stationary probability 
to find a random walker on a node of degree fc in a con¬ 
nected and undirected network with Mq links is given by 


k 

~ 2M^’ 


( 1 ) 


Using this relation, the probability hk^w) that w walkers 
are observed on a node of degree k is presented by 


hk{w) 





( 2 ) 


where Wq is the total number of walkers in the network. 
This leads a natural definition of the node capacity qk of 
a node of degree k as 


qk = {w)k + mak, 


( 3 ) 


where {w)k and Uk are the average and the standard devi¬ 
ation of the binomial distribution hk{w), which are given 
by {w)k = Wopk and ak = \/WoPfe(l -pk), respectively, 
and m is a real positive parameter which characterizes 
the node tolerance to load. Since the overload probabil¬ 
ity Fwq (k) of a node of degree k is the probability of w 
to exceed qk, F\^g{k) is given by summing up the distri¬ 
bution function Eq. m over w larger than qk. Thus, we 
have [2l| 


Fwoik)= E 

'W=lqk\+1 ^ 2 

= Ik/ 2 Mo{[Qk\ + ^,Wq — [qk\), (4) 

where Ip{a, b) is the regularized incomplete beta function 
[ 2 ^ and \_x\ is the greatest integer not greater than x. It 
is important to pay attention to the fact that the overload 
probability is a decreasing function of degree k |21|. 

Based on the above overload probability, we model the 
cascade process of overload failures as follows: 

(i) Prepare an initial connected, uncorrelated, and 
undirected network Qq with Nq nodes and Mg links, in 


(5) 

where Mr is the total number of links in the network Qr 
and r is a real positive parameter. 

(iii) Calculate the overload probability of every node, 
and remove nodes from Qr with this probability. 

(iv) Repeat (ii) and (iii) until no node is removed in 
the procedure (iii). 

In the procedure (ii), the total load Wr is reduced in ac¬ 
cordance with the reduction of the network size. In actual 
cases of cascading failures, the total load is often reduced 
to some extent during a cascade process to prevent the 
breakdown of the network function. For instance, when 
a problem arises in a power supply system due to a natu¬ 
ral disaster, the temporary restriction of electricity use is 
sometimes introduced to avoid a large-scale breakdown of 
the power grid, as in the case of the Tohoku earthquake 
and tsunami in 2011 (T7j |. In addition, when a company 
goes bankrupt, a large-scale chain bankruptcy could be 
prevented by the reduction of the total debt (loads) on 
the transaction network by means of, for example, a spe¬ 
cial low-interest lending facilities for companies having 
business relationships with the bankrupt company or the 
injection of taxpayers’ money. The quantity Wr given 
by Eq. ([5]) represents such a reduction of the total load 
during the cascade process. The exponent r character¬ 
izes how quickly the total load decreases with decreasing 
the network size, which we call hereafter the load re¬ 
duction parameter. Although the initial network Qq is 
connected, the network Qr at cascade step r is not nec¬ 
essarily connected. For an unconnected network Qr, we 
assume that Wr random walkers are distributed to com¬ 
ponents in proportion to the numbers of links in these 
components. Namely, the load allocated to the a-th com¬ 
ponent is given by 




( 6 ) 


where is the number of links in the a-th component. 

In the procedure (iii), the overload probability dur¬ 
ing the cascade process cannot be calculated straightfor¬ 
wardly by Eq. ®. First, the degree fc of a node in the 
network Qr at step r is not the same as its initial degree 
fcg. Since the probability to find a random walker on a 
node in Qr is a function of the present degree k of the 
node while the node capacity is determined by its initial 
degree fcg, the overload probability is presented by Eq. ® 
with replacing qk by qkg ■ Moreover, when the network Qr 
is unconnected, the overload probability must be calcu¬ 
lated in each component, because random walkers in a 
component cannot move beyond the component. There¬ 
fore, the overload probability of a node of degree fc, whose 
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initial degree is fcg, in the a-th component is given by 

Fw^iko,k) = /fc/2M“ (L9feo(W^o)J + 1 ,MKt “ l<lko0^o)\), 

(7) 

where M“ and are again the total number of links 
and the number of random walkers in the a-th component 
of Qt, respectively. 

The robustness of a network against cascading overload 
failures described above is evaluated by the relative size of 
the giant component S[ = N{/No, where N{ is the number 
of nodes in the largest component of the network Q{ at 
the final stage of the cascade process. More specifically, 
the robustness of the network is measured by the load 
reduction parameter Tc above which the relative size S{ 
becomes finite. A network providing a smaller Tc can 
be regarded to be more robust in the sense that there 
exists a giant component even if the total load is slowly 
reduced in accordance with the reduction of the network 
size during the cascade process. 

III. SIZE OF THE GIANT COMPONENT 

In order to assess the robustness of a network, we cal¬ 
culate the relative size St of the giant component in the 
network Qt at cascade step r. The calculation of Sr 
requires information on the load IT“ to obtain the over¬ 
load probability Fw^iko,k) for each component. It is, 
however, difficult to find an analytical expression of IT“. 
Thus, we assume that random walkers can jump to other 
components with a small probability, which enables us to 
estimate the stationary probability to find a walker on a 
node of degree k by 


instead of pg = A:/2M“, and the overload probability by 


average {w}'/. of this distribution function coincides with 
the average {w)k of h{w;W^,p'^). Although the stan¬ 
dard deviation cr^ = (w)/c(l — fc/2MT-) of h{w] WrjPk) 

is larger than Ok of h(w] IT“,p^), cr^ is not very different 
from cTfe if M“ (thus Mt) is large enough. Therefore, both 
distribution functions h{w;WT,pf) and h{w;WT,Pk) 
with the same average and similar widths are close to 
each other. Because of this similarity, F\Y^{ko,k) = 
Y.w>lqkai+iHf^l^?^Pk) that leads Eq. © can be ap¬ 
proximated by FwT{ko,k) = E™>[9fcj+i tAr.Pfc) 
giving Eq. (jH). 

We calculate the relative size St of the giant compo¬ 
nent by using the generating function method [1^. To 
this end, here we slightly modify the rule of the procedure 
(ii) described in Sec. II so that overloaded nodes are not 
removed but left in the system as zero-degree nodes for 
which random walkers never visit. This does not influ¬ 
ence any results in this work, except that the total num¬ 
ber of nodes in Qt remains constant at Nq, which makes 
the theoretical treatment easier. In spite of this modi¬ 
fication, we will continue to use the expression “remove 
a node” for simplicity, but this actually means “remove 
all links from a node”. The generating function method 
then enables us to calculate St if the network Qt is un¬ 
correlated and the degree distribution function PT{k) of 
Qt is given. In order to estimate PT{k)^ we introduce the 
probability nT-(fco,/c) that a randomly chosen node has 
the degree k and the initial degree feg. Obviously, the 
probability nT-(A:o, k) is related to the degree distribution 
PT{k) through 

PT{k)=Y,^r{ko,k). (11) 

ko^k 

As a special case, we have, at r = 0, 

no(fco, k) = Po{k)dkok, (12) 


FwT{ko,k) = Ik/2MT {[qko{Wo)\ +1 ,Wt- [qko{Wo)\), 

(9) 

instead of Eq. ©. This simplification does not largely 
change the overload probability if the number of links in 
a component is large enough. The reason is the following. 
The probability that w walkers are found on a degree-A: 
node in the a-th component containing IT“ walkers is 
written as 


because k is always equal to ko in the initial network 
Qq. To obtain the probability IlT{ko, k), we further in¬ 
troduce the probability (j)T{k) that the overload failure 
occurs on a node connected to a node of degree k in the 
network Qt- Considering that n^(A:o, k)/PT{k) represents 
the probability of a degree-fc node chosen randomly from 
Qt to have the initial degree fco, the probability 4>T{k) is 
expressed as 


^ypini-pif^-'^, ( 10 ) 

if random walkers are confined in the component. The 
average and the standard deviation of this binomial dis¬ 
tribution function of w are {w)k = W^P^ = WTk/2MT 
and (Tfc = \J (ly)fc(l — A:/2M“), respectively, where we 
used Eq. ([6]) and the relation p^ = fc/2M“. In the case 
that random walkers are allowed to jump between com¬ 
ponents, on the other hand, the probability to find w 
walkers on a degree-fc node is given by h(w] WT,Pk)- The 


ko 

•Prik) = P-rik'lk) 

ko k' — l 


nr(fco,fc') 

Pr{k') 


FwT{ko,k'), 


(13) 


where PT{k'\k) is the conditional probability that a node 
of degree k is connected to a node of degree k’. Since 
the network Qq has no degree correlations, Qt generated 
by removing nodes randomly from Qq with a probabil- 
ity depending only on the degree is also uncorrelated 
[29|. Therefore, the conditional probability is presented 
by Pr(A:'|fc) = k'PT{k')/{k)T, where (fc),- is the average 






4 


degree of Qr, and (pr is written as equal to the probability that a randomly chosen node 

from Qr-i has the initial degree fcp and the degree of 
k') , this node becomes fc during the cascade from C/t —1 to C/t- 

— Fw^iko,k'), (14) we can set up the master equation for 11 ,-(/cq, k) as 

which is independent of k. The probability Hriko^k) is 

I 


ko 

= EE 

ko 




k'n^jko, 

{k)r 


n,-(fco,fc) = ^ Il-r-l 

k'>k 




<(>r-i)^[l - Fw^_i{ko,k')] + SkoFw^-i{ko,k') 


(15) 


The right-hand side of this equation represents the prob¬ 
ability that a degree-fc' node in Gt-i becomes a node of 
degree k. The first term describes the situation that the 
degree-fc' node survives and k' — k nodes adjacent to this 
node are removed by overload failures. The second term 
stands for the case that the degree-fc' node itself experi¬ 
ences an overload failure and becomes a zero-degree node. 
Solving numerically Eq. (fT5)l with the aid of Eq. da, 
we can obtain the time evolution of the the probability 
Ur{ko,k) and the degree distribution Pr{k) by Eq. (ITTI) . 
The relative size St of the giant component at cascade 
step T is calculated by [2^ 

= 1-^P.,(fc)^^^ (16) 

k 

where u is the smallest non-negative solution of the tran¬ 
scendental equation, 


u = 


(17) 


and is the generating function of the remaining 

degree distribution, which is defined by 


G''^\x) = —Y,{k + l)PT{kSl)x\ (18) 


According to the procedure (iv) stated in Sec. m the 
cascade process must be terminated when no node ex¬ 
periences overload failures. In an actual calculation, we 
stop the cascade process at the step r satisfying the con¬ 
dition 


'^Fw^{ko,k)IlT{ko,k) < (19) 

kto 

This condition implies that the expectation number of 
overloaded nodes becomes less than unity. 

It should be noted that the above formalism is based 
on the approximation that random walkers are allowed to 
jump over components with a small probability. In order 
to evaluate the accuracy of this approximation, we com¬ 
pare the time evolutions of St for an Erdos-Renyi random 
graph calculated by both the method explained in this 



FIG. 1. (Color online) Relative size St of the giant com¬ 
ponent as a function of the cascade step r. The initial net¬ 
work is an Erdos-Renyi random graph with 5, 000 nodes and 
12, 500 links, which implies that the initial average degree is 
{k)o = 5.0. The symbols indicate the results obtained by nu¬ 
merical simulations following the cascade process from (i) to 
(iv) described in Sec. El and averaging over 50 network real¬ 
izations. Error bars represent one standard deviation from 
the mean values. The lines show the results calculated by the 
analytical method explained in Sec. ED The total number of 
walkers at the initial cascade step is Wo = 25, 000. The val¬ 
ues of the node tolerance parameter m and the load reduction 
parameter r are displayed in the hgure. 


Section and numerical simulations following the cascade 
process from (i) to (iv) described in Sec. HD The results 
shown in Fig. [T] indicate that St calculated by Eq. (fTCl) 
agrees quite well with the simulation result, which reflects 
the high accuracy of the approximation. 


IV. RESULTS 

We assess the robustness of networks against cascad¬ 
ing overload failures by computing the relative size S{ 
of the giant component at the final stage of the cascade 
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FIG. 2 . (Color online) Relative size St of the giant compo¬ 
nent at the final stage of the cascade process as a function 
of the load reduction parameter r, for an ER random graph 
(black dashed line and filled squares) and an SF network (red 
solid line and filled circles). The SF network has the degree 
distribution given by Eq. ([ 20 } with 7 = 2 . 5 . The lines repre¬ 
sent the results calculated by the analytical method explained 
in Sec. imi The symbols indicate the simulation results ob¬ 
tained by following numerically the cascade process from (i) 
to (iv) described in Sec. m and averaging over 50 samples. 
Only error bars larger than the symbol size are shown. In 
this calculation, we set No = 10 "^, Mo = 2 x 10 "*, m = 2 . 0 , 
and Wo ~ 2Mo. 


process as a function of the load reduction parameter r. 
In this work, we examined two types of initial networks, 
namely, homogeneous Erdos-Renyi (ER) random graphs 
and scale-free (SE) networks with inhomogeneous degree 
distributions. An SF network is generated by the con¬ 
figuration model [ 1 ^ with the degree distribution given 
by 

Po(fc) = ( ^ (20) 

[O for fc < fcmin, 

where d is a real positive parameter and 7 is the exponent 
describing the asymptotic power-law form of the degree 
distribution, i.e., Po{k) oc k~'^ for k ^ d. In Eq. (1^ . 
^min is the minimum degree, and C is the normalization 
constant. The parameters d and fcmin can control the av¬ 
erage degree (fc)o for a fixed value of 7 . Throughout this 
paper the minimum degree is fixed at fcmin = 2. Figure 
shows the r dependence of St for an ER random graph 
(black dashed line and filled squares) and for an SF net¬ 
work (red solid line and filled circles). Both networks 
have Nq = 10"'’ nodes and Mq = 2 x 10"' links, which 
implies that the initial average degree is (fc)o = 4.0. The 
initial total load is chosen as Wq = 2Mq and the node 
capacity qk is determined by Eq. (jS]) with m = 2.0. For 
the SF network, the exponent 7 is set as 7 = 2.5 and d 


FIG. 3. (Color online) Scaling plot of S'f(r, No) for ER ran¬ 
dom graphs of several sizes {No = 2"", 2'®, 2'®, and 2®®). The 
average degree of these networks is set to be (fc)o = 4.0. The 
calculations have been done under the condition of Wo = 
No{k)o and m = 2.0. In this plot, we choose rc = 1.84, 
P = 1.27, and l/v* = 19.45. Two branches correspond to the 
percolating and non-percolating phases. 


is tuned to satisfy (fc)o = 4.0. The solid and dashed lines 
in Fig. [5] indicate the results calculated by the analytical 
method explained in Sec. uni whereas the symbols rep¬ 
resent Sf obtained by simulating numerically the cascade 
process from (i) to (iv) described in Sec. m For both 
types of networks, the analytical results agree quite well 
with the numerical ones. As shown in Fig. [51 there ex¬ 
ists a value rdNo) below which 5'f = 0 and above which 
5'f > 0. This implies that a global cascade of overload 
failures occurs if the total load Wr is reduced during the 
cascade more slowly than the reduction scheme Eq. 
with r = rdNo) while a finite fraction of nodes survives 
the cascade if the reduction of Wr is faster. The fact 
that the value of rc{No) for the SF network is smaller 
than that for the ER random graph suggests that SF 
networks are more robust to cascading overload failures 
than ER random graphs. This tendenw is opposite to 
what was predicted by previous works |9l-[l^. 

The above conclusion is, however, obtained for finite- 
size networks. It is crucial to study the robustness of 
infinitely large SF networks against cascading overload 
failures. If the change of St at r = rdNo) provides a crit¬ 
ical transition when Nq goes to infinity, we can use the 
finite-size scaling analysis to calculate Tc in the thermo¬ 
dynamic limit. According to the finite-size scaling theory, 
the relative size of the giant component in the network 
Qf at the final stage of the cascade process starting from 
the initial network Qo of size Nq is expressed as 

Stir, No) = k - ^cl) , (21) 

where the correlation volume exponent ly* characterizes 
the divergence of the number of nodes within the 
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FIG. 4. (Color online) Nb{Ib) for the giant components in 
C/f’s obtained numerically by simulating the cascade process 
with r = rc{No) starting from an ER random graph (black 
filled squares) and an SF network with 7 = 3.5 (red filled 
circles). The number of nodes and the average degree of both 
of the initial networks are No = 10“' and (fc)o = 4.0, respec¬ 
tively. The calculations have been done under the condition 
of Wo = No(k)o and m = 2.0. The results are averaged over 
100 samples. The longitudinal axis is rescaled by Nb{Ib) at 
/b = 1. The dashed lines are merely guides to the eye. 


FIG. 5. Critical load reduction parameter as a function of 
the exponent 7 in the thermodynamic limit. In this calcula¬ 
tion for SF networks, we tune the parameter d in Eq. (I 20 I 1 to 
satisfy (fc)o = 4.0 for various values of 7 . The values of are 
estimated through the finite-size scaling analysis for S'f(r, No) 
calculated by the method described in Sec. [ITT] under the con¬ 
dition of Wo = {k)oNo and m = 2.0. The horizontal dashed 
line is the result for the ER random graph (7 —>■ 00 ) with 
{k)o = 4.0. 


correlation volume as oc \r — rc\~'^ in the infinite 
system, /? is the critical exponent for the approach to zero 
of Si{r, 00 ), Tc is the critical load reduction parameter in 
the thermodynamic limit, and $(a:) is a scaling function. 
Therefore, if the suitable values of the parameters v*, 
13, and Tc are selected, the quantity S{{r,No) as a 

function of |?' — rd collapses onto a single curve for 

various values of r and Nq. Figure [3] shows such a plot 
for ER random graphs of different sizes with the use of 
the r dependence of Sf{r,No) calculated by the method 
explained in Sec. Ell and the best-fit values of v*, (3, and 
Tc. The fact that all data collapse onto a single curve 
implies that the transition aX r = can be considered 
as a critical phenomenon. Similar scaling behaviors have 
been confirmed for SF networks. 

The criticality of the network C/f at r = Tc has also been 
conhrmed by the fractal property of the giant component 
in Qf. The fractality of complex networks is widely in¬ 
vestigated by the box covering algorithm [El - ls^ . If the 
minimum number Nb(Ib) of subgraphs of radius Ib re¬ 
quired to cover a given connected network satisfies the 
relation 

( 22 ) 

the network is considered to be fractal with the fractal 
dimension de [3l|. We calculated, by using the compact- 
box-burning algorithm A^b(^b) for giant components 


included in networks Qi at the final stage of the cascade 
process with r = Tc^Nq) starting from both an ER ran¬ 
dom graph and SF network. The giant components are 
obtained numerically by simulating the cascade process 
described in Sec. El The results shown in Fig. |4] indi¬ 
cate that the structures of these giant components ex¬ 
hibit the fractal nature, which supports the criticality of 
C/f at r = Tc. The fractal dimension for the ER random 
graph is de = 1-54 ± 0.01, while de = 2.13 ± 0.02 for 
the SF network with 7 = 3.5. These values of de are 
different from those of giant components at the critical 
point of the conventional percolation with random node 
removals, which are de = 2 for ER random graphs (or 
SF networks with 7 > 4) and de = (7 — 2)/(7 — 3) for 
SF networks with 3 < 7 < 4 [s^. Such a discrepancy is, 
of course, due to the difference in ways of node removals. 
In the cascade process, nodes that will be removed at 
the cascade step t depends strongly on nodes removed 
at r — 1 , as in the case of a fire spread for which a por¬ 
tion remaining unburned has a different structure from 
that of a survival from random removals. 

The critical load reduction parameter Tc calculated by 
the hnite-size scaling method is plotted in Fig. [5] as a 
function of the exponent 7 of SF networks. In this hg- 
ure, fc for the ER random graph is also indicated by the 
horizontal dashed line. The results only for 7 > 5 are 
presented here, because for a small value of 7 the com¬ 
putation of Ilr{ko,k) requires a long CPU time due to 
the increase of the maximum degree fc^ax of the initial 
network Qq associated with a decrease of 7 as is given by 
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fcmax oc -/Vg^^. We see from Fig.[5]that Tc is an increasing 
function of 7 . This implies that the enhancement of the 
SF property by decreasing 7 makes networks robust to 
cascading overload failures even in the thermodynamic 
limit. It should be emphasized that rc is always posi¬ 
tive for any 7 . This is because in any network a cas¬ 
cade of overload failures under r = 0 never stops until 
the network collapses c omp letely. Our result contrary to 
previous predictions [9I41 It| comes from the fact that the 
overload probability is a decreasing function of degree 
k. At the first step of the cascade process, nodes with 
small degrees are more likely to be removed according 
to the overload probability Fy/^{k). A low-degree node 
tends to be connected to a node with large degree in 
an uncorrelated SF network. The overload probability 
of a high-degree node adjacent to the low-degree node 
which was removed at the first step becomes smaller at 
the second cascade step than its initial overload proba¬ 
bility. Therefore, nodes with relatively small degrees are 
again preferentially removed also at the second step, and 
so on. It is obvious that SF networks are robust against 
preferential removals of low-degree nodes . 

Our model differs from previous models of cascading 
overload failures in two points. One is that overload fail¬ 
ures in our model occur when fluctuating loads exceed 
the capacities predetermined for nodes, while failures are 
caused by non-fluctuating loads (or average values of fluc¬ 
tuating loads) exceeding the node capacities in previous 
models [in. The other difference is in the dynamics of 
loads on a network. In many of previous works , the 

load is defined by flow passing through the shortest path 
between a pair of nodes. This type of loads describes, for 
example, a packet flow in the Internet and a traffic flow 
in a road system. On the other hand, loads in our model 
move randomly on a network, as in the case of a flow 
of debt in a corporate transaction network. It is impor¬ 
tant to clarify whether the robustness of SF networks is 
caused by the fluctuation of load or by its random walk 
behavior. 

In order to identify the origin of the robustness, 
namely, the origin of Fwo{k) being a decreasing function 
of degree k, we calculate the overload probability for fluc¬ 
tuating loads imposed by a shortest-path flow. Fluctu¬ 
ating loads carried b y fl ow along the shortest paths have 
been argued in Ref. [1^ to explain the relation between 
the average flux and the fluctuations. According to this 
work, we newly define the load Wi of the node i by 


Wi = 


E 


(») 


(23) 


where Vu is a set of u node pairs that are randomly se¬ 
lected from {Nq — 1) (iVo — 2)/2 node pairs in Qq excluding 
the node i, ajji is the total number of shortest paths be¬ 
tween the pair (j,j') in 14, and is the number of 

those paths that pass through i. Since the quantity rep¬ 
resented by the right-hand side of Eq. (USD is equivalent 
to the betweenness centrality bi if I 4 is chosen as the set 



FIG. 6 . Degree dependence of the overload probability Fu (k) 
for several values of the node tolerance parameter m. The 
lines from top to bottom represent the results for m = 2 . 0 , 
3.0, 4.0, 5.0, and 6.0, respectively. The size of the set Vu 
is fixed at m = 500. The load distribution hk{w) required 
to calculate Fu{k) is numerically evaluated by computing the 
partial betweenness centrality for each of 10 “^ different sets Vu 
selected from each of 30 SF networks with 10"^ nodes. The 
degree distribution of these SF networks is given by Eq. (Uni) 
with 7 = 2.5 and d = 0.152 providing (k) = 4.0. 


of all node pairs, namely, u = {Nq — l)(A^o ~ 2 )/ 2 , we 
call the above quantity Wi the partial betweenness cen¬ 
trality. The partial betweenness centrality of a node de¬ 
pends on which node pairs are selected for the set I 4 , 
and thus the load Wi fluctuates in accordance with the 
choice of 14. For a large number of different sets 14 with 
1 <C u <C (IVo —l)(-^o — 2 )/ 2 , the distribution hi{w) of the 
load can be calculated for each node. As the conventional 
betweenness centrality bi is strongly correlated to the de¬ 
gree ki in an uncorrelated network [dll, jd^ , the partial 
betweenness centrality of a node is also expected to have 
a correlation with the degree of the node. In fact, we 
have confirmed numerically that two distribution func¬ 
tions hi{w) and hj{w) of load on the nodes i and j which 
have the same degree are close to each other. There¬ 
fore, we can define the degree-dependent load distribu¬ 
tion hk{w) for the partial betweenness centrality model, 
which corresponds to the load distribution hk{w) given 
by Eq. dH) for the random walker model. The node capac¬ 
ity qk is also defined by Eq. m with the average load {w)k 
and the standard deviation Ufe calculated by the distribu¬ 
tion function hk{w). Finally, as in the case of the random 
walker model, the overload probability is defined as the 
probability that the load of a node of degree k exceeds its 
capacity qk, which is calculated by hk(w)dw. Since 

hk(w) depends on u, we denote this overload probability 
by Fu{k) instead of Fwo{k). 

Figure [S] shows the degree dependence of F„(fc) for sev¬ 
eral values of the node tolerance parameter m. We see 









TABLE 1. Robustness/fragility of SF networks against cas¬ 
cading overload failures induced by several types of loads. 
The fragility for shortest path flow without fluctuation is a 
consequence of the previous work [^. 


Load flow 

Without fluctuation 

With fluctuation 

Shortest path 

Fragile 

Robust 

Random walk 

Fragile 

Robust 


that the overload probability F„(fc) is a decreasing func¬ 
tion of k for any m. This implies that SF networks are 
robust also against cascading overload failures induced 
by fluctuating shortest-path flow. It is not surprising 
that the degree dependences of Fyyg{k) and Fu{k) show 
a similar tendency, because it has been reported that the 
couplings between the fluctuations and the average of the 
number of random walkers and the partial betweenness 
centrality on individual nodes obeys the same scaling law 
[T^ . We can also consider a situation that cascading 
failures are caused by non-fluctuating loads by random 
walkers. In this case, an overload failure occurs when the 
average number of walkers {w)k = Wok/2Mo on a node 
of degree k exceeds the node capacity depending on {w)k, 
where Wq is the total number of walkers and Mq is the 
number of links in the initial network. Since {'w)k is pro¬ 
portional to k and the degree of a node correlates closely 
with the betweenness centrality of the node [4l|, , the 

property of cascading failures induced by non-fluctuating 
random walking loads is essentially the same as that by 
loads of the betweenness centrality 0. Therefore, SF 
networks are fragile to such cascading failures. All the 
above arguments can be summarized as shown in TablelH 
From this Table, we can conclude that the robustness of 
SF networks in our model is a consequence of the prop¬ 
erty that failures are caused by extreme values of fluc¬ 
tuating loads, but not concerned with the random walk 
behavior of loads. 


V. CONCLUSIONS 

We have studied the robustness of scale-free (SF) net¬ 
works against cascading overload failures induced by ex¬ 
treme values of fluctuating loads that exceed the node 
capacities. In our model, temporally fluctuating loads 
are treated as random walkers on a network, for which 
the stationary overload probability has been studied by 
Kishore et al. [2l|. At the first stage of the cascade, nodes 
are removed from the initial network with this overload 
probability, and the redistribution of loads in the dam¬ 
aged network causes iteratively subsequent failures ac¬ 
cording to the updated overload probability until no node 
is expected to be removed. During the cascade process, 
the total load is reduced in response to the decrease of the 
number of links in the network under the cascade. How 


quickly the total load is reduced is characterized by the 
load reduction parameter r. The robustness of a network 
is measured by the critical load reduction parameter 
above which the relative size Si of the giant component 
at the final cascade stage is finite. We present a formu¬ 
lation to calculate S{ by using the master equation for 
the probability n,-(A:o, k) of a node in the network Qt at 
cascade step r to have the present degree k and the ini¬ 
tial degree ko and by applying the generating function 
method. Our results for SF networks with degree dis¬ 
tribution Po{k) ^ k~^ at large k show that Cc increases 
with the exponent 7 , which implies that SF networks are 
robust against cascading overload failures in our model 
as opposed to previous works [iSl. The robustness of 
SF networks is explained by the property of the overload 
probability of being a decreasing function of the degree, 
which is a consequence of the load fluctuations but not 
concerned with the random walk behavior of loads. 

In our model, the total load does not fluctuate through¬ 
out the cascade process though the local load on a node 
fluctuates. If the total load also fluctuates temporally 
and the time scale of the fluctuation is faster than that 
of the cascade process, the overload probability must be 
different from Eq. o or ®. However, we can expect 
that the overload probability remains to be a decreas¬ 
ing function of the degree if the magnitude of the total 
load fluctuation A IF is much smaller than the average 
total load (IF), which guarantees the robustness of SF 
networks to cascading overload failures. This is because 
the standard deviation ak of the fluctuating load w on a 
node of degree k is proportional to the square root of the 
average load {w)k if AIF <C (IF), as in the case of the 
present work. It has, however, been reported that ak be- 
comes proportional to {w)k when AIF/(IF) approaches 
to 1 [ 1 ^. The change in the property of the load fluc¬ 
tuations may alter drastically the degree dependence of 
the overload probability and hence the robustness of net¬ 
works. It is thus important to study how robust SF net¬ 
works are against cascading failures under large flnctua- 
tions of the total load. Furthermore, although we present 
the formulation to calculate the relative size of the giant 
component Si in this work, an analytical expression for 
the critical load reduction parameter Tc is not found. In 
addition to solving this problem, the identification of the 
universality class of the present model is also a subject 
for future work. 


ACKNOWLEDGMENTS 

This work was supported by a Grant-in-Aid for Sci¬ 
entific Research (Nos. 25390113 and 14J01323) from the 
Japan Society for the Promotion of Science. Numerical 
calculations in this work were performed in part on the 
facilities of the Supercomputer Center, Institute for Solid 
State Physics, University of Tokyo. 







9 


[1] D. Krol, D. Fay, and B. Gabrys, Propagation Phenomena 
in Real World Networks (Springer, Heidelberg, 2015). 

[2] L. da F. Costa, O. N. Oliveira, G. Travieso, F. A. Ro¬ 
drigues, P. R. V. Boas, L. Antiqueira, M. P. Viana, and 
L. E. C. da Rocha, Adv. Phys. 60, 329 (2011). 

[3] R. Cohen, K. Erez, D. ben-Avraham, and S. Havlin, 
Phys. Rev. Lett. 85, 4626 (2000). 

[4] R. Cohen, K. Erez, D. ben-Avraham, and S. Havlin, 
Phys. Rev. Lett. 86 , 3682 (2001). 

[5] M. E. J. Newman, Phys. Rev. Lett. 89, 208701 (2002). 

[ 6 ] P. Holme, B. J. Kim, C. N. Yoon, and S. K. Han, Phys. 
Rev. E 65, 056109 (2002). 

[7] S. Boccaletti, V. Latora, Y. Moreno, M. Chavez, and D.- 
U. Hwang, Phys. Rep. 424, 175 (2006). 

[ 8 ] S. N. Dorogovtsev, A. V. Goltsev, and J. F. F. Mendes, 
Rev. Mod. Phys. 80, 1275 (2008). 

[9] A. E. Motter and Y.-C. Lai, Phys. Rev. E 66 , 065102(R) 

( 2002 ). 

[10] P. Holme, Phys. Rev. E 66 , 036119 (2002). 

[11] P. Crucitti, V. Latora, and M. Marchiori, Phys. Rev. E 
69, 045104(R) (2004). 

[12] J. J. Wu, H. J. Sun, and Z. Y. Gao, Physica A 386, 407 
(2007). 

[13] Z. J. Bao, Y. J. Cao, L. J. Ding, and G. Z. Wang, Physica 
A 388, 4491 (2009). 

[14] Y. Xia and D. J. Hill, Physica A 389, 1281 (2010). 

[15] B.-L. Don, X.-G. Wang, and S.-Y. Zhang, Physica A 389, 
2310 (2010). 

[16] M. A. de Menezes and A.-L. Barabasi, Phys. Rev. Lett. 
92, 028701 (2004). 

[17] S. Meloni, J. Gomez-Gardenes, V. Latora, and 
Y. Moreno, Phys. Rev. Lett. 100, 208701 (2008). 

[18] Z. Zhou, Z.-G. Huang, L. Huang, Y.-C. Lai, L. Yang, and 
D.-S. Xue, Phys. Rev. E 87, 012808 (2013). 

[19] P. Embrechts, C. Klueppelberg, and T. Mikosch, Mod¬ 
eling Extremal Events: For Insuranee and Finance 
(Springer, Berlin, 1997). 

[20] Fluctuation-induced cascading failures have been dis¬ 
cussed in D. Heide, M. Schafer, and M. Greiner, Phys. 
Rev. E 77, 056103 (2008). 

[21] V. Kishore, M. S. Santhanam, and R. E. Amritkar, Phys. 


Rev. Lett. 106, 188701 (2011). 

[22] V. Kishore, M. S. Santhanam, and R. E. Amritkar, Phys. 
Rev. E 85, 056120 (2012). 

[23] S. Mizutaka and K. Yakubo, Phys. Rev. E 88, 012803 
(2013). 

[24] S. Mizutaka and K. Yakubo, J. Complex Networks 2, 413 
(2014). 

[25] J. D. Noh and H. Rieger, Phys. Rev. Lett. 92, 118701 
(2004). 

[26] M. Abramowitz and 1. A. Stegun, Handbook of Mathe¬ 
matical Functions (Dover, New York, 1964). 

[27] H. Arikawa, Y. Cao, and S. Matsumoto, Energy Res. Soc. 
Sci. 2, 12 (2014). 

[28] M. E. J. Newman, S. H. Strogatz, and D. J. Watts, Phys. 
Rev. E 64, 026118 (2001). 

[29] A. Srivastava, B. Mitra, N. Ganguly, and F. Peruani, 
Phys. Rev. E 86, 036106 (2012). 

[30] M. Molloy and B. Reed, Random Struct. Algor. 6, 161 
(1995). 

[31] C. Song, S. Havlin, and H. A. Makse, Nature (London), 
433, 392 (2005). 

[32] C. Song, S. Havlin, and H. A. Makse, Nat. Phys. 2, 275 
(2006). 

[33] K.-I. Goh, G. Salvi, B. Kahng, and D. Kim, Phys. Rev. 
Lett. 96, 018701 (2006). 

[34] J. S. Kim, K.-I. Goh, G. Salvi, E. Oh, B. Kahng, and 
D. Kim, Phys. Rev. E 75, 016110 (2007). 

[35] F. Kawasaki and K. Yakubo, Phys. Rev. E 82, 036113 

( 2010 ). 

[36] S. Furuya and K. Yakubo, Phys. Rev. E 84, 036118 

( 2011 ). 

[37] Y. Sun and Y. Zhao, Phys. Rev. E 89, 042809 (2014). 

[38] C. Song, L. K. Gallos, S. Havlin, and H. A. Makse, 
J. Stat. Mech.: Theory Exp. (2007), P03006. 

[39] R. Cohen and S. Havlin, Physica A 336, 6 (2004). 

[40] A. A. Moreira, J. S. Andrade, Jr., H. J. Herrmann, and 
J. O. Indekeu, Phys. Rev. Lett. 102, 018701 (2009). 

[41] K.-I. Goh, B. Kahng, and D. Kim, Phys. Rev. Lett. 87, 
278701 (2001). 

[42] M. Barthelemy, Eur. Phys. J. B 38, 163 (2004). 



